Neutrino oscillations with three flavors in matter of varying density 



O 

o 

(N 



Tommy Ohlsson* 

Institut fur Theoretische Physik, Physik Department, Technische Universitat Miinchen, 
James-Franck-Strajie, DE-85748 Garching bei Miinchen, Germany 

and 

Division of Mathematical Physics, Theoretical Physics, Department of Physics, Royal Institute of Technology, 

SE-100 44 Stockholm, Sweden 

Hakan Sncllman^ 

Division of Mathematical Physics, Theoretical Physics, Department of Physics, Royal Institute of Technology, 

SE-100 44 Stockholm, Sweden 
(February 1, 2008) 

In this paper, we discuss the evolution operator and the transition probabilities expressed as 
functions of the vacuum mass squared differences, the vacuum mixing angles, and the matter density 
parameter for three flavor neutrino oscillations in matter of varying density in the plane wave 
approximation. The applications of this to neutrino oscillations in a model of the Earth's matter 
density profile, step function matter density profiles, constant matter density profiles, linear matter 
density profiles, and finally in a model of the Sun's matter density profile are discussed. We show that 
for matter density profiles, which do not fluctuate too much, the total evolution operator consisting 
of n operators can be replaced by one single evolution operator in the semi-classical approximation. 
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In previous papers we have given analytic expressions for the three flavor neutrino oscillation evolution 

operator and the transition probabilities in presence of constant matter densities expressed in the vacuum mixing 
■ matrix elements and the neutrino energies or masses, i.e., incorporating the so-called Mikheyev-Smirnov-Wolfenstein 
(MSW) effect Here we will discuss the application of this to realistic matter density variations in a "semi- 

classical" approximation based on our previous results. This allows a simple and efficient calculation of neutrino 
oscillations in media of varying densities. We compare this approximate formula with a numerical simulation in a 
multi-step model. We will as before assume that the CP phase S is equal to zero. Thus, the neutrino mixing matrix 
is real. The semi-classical approximation for three neutrino flavors, we believe, is a unique part of our investigation. 

Previous work on models for three flavor neutrino oscillations in matter for constant matter density includes works 
of Barger et al. || , Kim and Sze , and Zaglauer and Schwarzer JtJ . Approximate solutions for three flavor neutrino 
oscillations in matter have been presented by Kuo and Pantaleone B and Joshipura and Murthy 0. Approximate 
treatments have also been done by Toshev and Petcov [ fTo|| , D'Olivo and Oteo have made contributions by using an 
approximative Magnus expansion for the time evolution operator pd|. Extensive numerical investigations for matter 
enhanced three neutrino oscillations have been made by Fogli et al. ||12|| . Studies of neutrino oscillations in Earth has 



I. INTRODUCTION 



been performed by several authors [pL3| |l8 1 



Neutrino oscillations for matter with linearly varying density have been treated by Petcov |19[ and Lehmann et al. 
pof . Osland and Wu [M have also solved the case for exponentially varying density. Matter enhanced two flavor 
neutrino oscillations with an arbitrary monotonic matter density profile have been studied by Balantekin and Beacom 
p2| using a uniform semi-classical approximation. See also Fishbane et al. p3fl for two flavor neutrino oscillations in 
matter of varying density. 
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II. THE EVOLUTION OPERATOR IN PRESENCE OF MATTER 



Let the flavor state basis and mass eigenstate basis be denoted by Bf = {\v a )}a=e,fj,,T and B m = {l^a)}^].; 
respectively. Then, the flavor states \v a ) € Bf can be obtained as a superpositions of the mass eigenstates \v a ) & Bm, 
or vice versa. The bases Bf and B m are of course just two different representations of the same Hilbert space Tt. 

In the present analysis, we will use the plane wave approximation to describe neutrino oscillations. In this approx- 
imation, a neutrino flavor state \v a ) is a linear combination of neutrino mass eigenstates |^ a )'s such that p4|| 

3 

K) = Y. U *M, (i) 

a=l 

where a = e, fi, r. In what follows, we will use the short-hand notations \a) = \v a ) and \a) = \v a ) for the flavor states 
and the mass eigenstates, respectively. 

The components of a state ip in flavor basis and mass basis, respectively, are related to each other by 



0/ = U^ m , (2) 



where 



1pe\ I tpl 

ip f ee (ip a ) ee | tpn e Bf and tp m = (ip a ) = ip 2 | € B m . 

4>r J \ 03 

A convenient parameterization for U = U (81,62,93) is given by ||{| 

C2C3 S3C2 
U = I —S3C1 — S1S2C3 C1C3 — S1S2S3 S1C2 ] , (3) 
S'i'S'3 — S2C1C3 —S1C3 — S2S3C1 

where Si = sin 9i and Ci = cos 6i for i = 1,2,3. This is the standard representation of the neutrino mixing matrix. 
The quantities 9i, where i = 1, 2, 3, are the vacuum mixing angles. Since we have put the CP phase equal to zero in 
the neutrino mixing matrix, this means that U* a = U aa for a = e, //, r and a = 1,2,3. 

In mass basis, the Hamiltonian ffl 1 for the propagation of the neutrinos in vacuum is diagonal and given by 

Ei 

H m = I E 2 1 , (4) 




E; 
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where E a = m 2 a + p 2 , a = 1, 2, 3, are the energies of the neutrino mass eigenstates |a), a = 1, 2, 3 with masses m a , 
a = 1, 2, 3. We will assume the three-momentum p to be the same for all mass eigenstates. 

When neutrinos propagate in ordinary matter, there is an additional term in the Hamiltonian J4? coming from the 
presence of electrons in matter B . This term, the potential term, is diagonal in flavor basis and is given by 




V f = A I ) EE ,I/V ; . (.'.) 

where 

A ee A(r) = ±V2G F NJr) ~ ±^=G F — p{r) 

V2 m N 

is the matter density parameter and Kf is the projector in flavor basis on the electron neutrinos. Here Gf is the 
Fermi weak coupling constant, N e is the electron density, win is the nucleon mass, and p is the matter density. The 
sign of the matter density parameter depends on weather we deal with neutrinos (+) or antineutrinos (— ). In mass 
basis, this piece of the Hamiltonian is V m = U~ 1 VfU, where U is again the neutrino mixing matrix. 

In the case when the neutrinos propagate through matter, as here, the Hamiltonian is not diagonal in either the mass 
basis or the flavor basis, and we have to calculate the evolution operator Uf(t) or Uf(L) = e - l -^ L = u e -^mLjj-i 
if we set t = L (L is the traveling (propagation) path length of the neutrinos.). 
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To do so it is convenient to introduce the traceless real symmetric matrix T defined by T = Jf m — (tr Jf m )I/3. The 
trace of the Hamiltonian in mass basis Jif m = H m + U~ 1 VfU is tr Jif m = E\ + E% + E 3 + A, and the matrix T can 



then be written as 



AUl x -hA+h (E12 + E 13 ) AU el U e2 AU el U t 



c:3 



T = (T ab ) = ( AU el U e2 AU 2 2 - \A + § (E 21 + E 23 ) AU e2 U e3 | , (6) 

AU el U e3 AU e2 U e3 AU 2 e3 -\A+\ (E 31 + E 32 ) 

where E ab = E a — E b . Of the six antisymmetric quantities E ab , where a,b = 1,2,3 and a ^ b, only two are linearly 
independent, since the E ab 's fulfill the relations E ba = -E ab and E 12 + E 23 + E 31 = oj] This means that the evolution 
operator in mass basis can be written as 



U m (L) = e~^ L = ^e~ lLT = tJ2e-* LX * —J— [(X 2 a + Cl )I + \ a T + T 2 ] , (7) 

where <j> = e ~ lL ( tr ^"i) / / 3 j \ a , a = 1, 2, 3, are the eigenvalues of the matrix T, and / is the 3x3 identity matrix.^ The 
coefficients cq, ci, and c 2 are all real and the eigenvalues A Q , a = 1, 2, 3, can be expressed in closed form in terms of 
these 

The evolution operator for the neutrinos in flavor basis is thus given by 

3 1 r i 

U^L^e-i^fL = Ue -iM' m L u -i = ^ e - iL x a (X 2 a + Cl )I + X a f + f 2 , (8) 

a=l + Cl L J 

where T = UTU^ 1 . Equation (^) is our final expression for Uf(L). 

Since Jtff = UJ^JJ , it is clear that T = J%f — (tr J0f)I /3 — J0f — (tr J^f m )I /3 due to the invariance of the 
trace under transformation of U. In fact, the characteristic equation is also invariant under transformation of U and 
therefore so are the coefficients Cq, c\, c 2 , and the eigenvalues Ai, A 2 , A3. However, the expression for Jiff is much 
more complicated than that for Jif m , which is the reason why we work with Jif m instead of Jiff. 

The formula (||) expresses the time (or L) evolution directly in terms of the mass squared differences and the 
vacuum mixing angles without introducing any auxiliary matter mixing angles. By dividing the density variation in 
small, approximately constant segments, and using this formula repeatedly in each segment, we can numerically study 
neutrino oscillations in matter with varying density. We will use this method as a standard test for the semi-classical 
approximation to the evolution operator that we study below. 



III. THE SEMI-CLASSICAL APPROXIMATION 

The formula given in Eq. Q is the evolution operator in flavor basis for constant matter density. To handle the 
case of varying matter density, let us divide the distance L from the source to the detector into N equidistant parts 
and put the index k on the eigenvalues A^, where a = 1, 2, 3 denote the three mass eigenstates. For any matter density 
profile p(r) we first make this profile discrete and introduce pu as the matter density in the interval rk—i < r < 
where k varies from 1 to N with tq — and — L. The length of each segment is then Ar^ = r& — r^-i — L/N. 
The evolution operator in mass basis from to L can then be written as the ordered product 

U m (L) = U m {r N ~ r N -i)U m (r N ^i ~ r N - 2 ) ■ ■ ■ U m (r 2 - ri)U m (rx - r ). (9) 

Note that the order of the U m (rk — rfc_i)'s are important, since these operators do not in general commute. 



1 Later, we will use the usual (vacuum) mass squared differences Am^ and Am^, instead of £21 and E32, which are related 
to each other by Am|i = 2E V E2\ and A7TI32 — 2E„E32, where E v is the neutrino energy. 

2 Using Cayley-Hamilton's theorem, the exponential of a matrix M, e M = X^S^Lo ~k\M n (infinite series), can be written as 
e M _ a n M n (finite series, iV is the dimension of M), where a„ (n = 0, 1, . . . , N — 1) are some coefficients to be 

determined. In this case, since N = 3 (the dimension of T is three), this means that we have e _l = aol — ia\ LT -a 2 L 2 T 2 , 
i.e., there are no higher power terms of T in e~ lLT than that of order two. 
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When N is large, each step is small and the exponent in the evolution operator U m (rk — r k -x) = U m (Ark) is small. 
We can then approximate this operator with 

U m (Ar k ) ~ e" iA ^, (10) 
where J^ns = H rn + AkK m and A k cx p/.. Inserting this into Eq. (^) gives 

Since the ^n's do not commute, the higher order terms have to be calculated with the time-ordering (here rather 
r-ordcring) operator. However, here wc will at first be satisfied with the lowest order result, which we call the semi- 
classical approximation. In this approximation, we retain only the terms proportional to Ar k = L/N, and thus, 
neglect the noncommutativity of the ^^'s for different k's. We can thus write 

U m {L) ~ e-^i (12) 
In the limit N — > oo, this gives the integral formula 

U m (L) = e^Io^^dr = 0(£) e -i/ o I 'T(r)dr ) (13) 

where T(r) is the traceless part of the Hamiltonian corresponding to the electron density at position r between and 
L and <f>(L) is the phase factor coming from the trace. 

For further discussion it is often convenient to retain the original form of the Hamiltonian and to use J4? m = 
H m + AU~ 1 K f U rather than T. Thus, when A = A(r), we obtain 

L 

J%n{r)dr = L(H m + A{L)K m ), (14) 



where 

A(L) = - [ A(r)dr 
L Jo 

is the average matter density along the baseline L and 

K m ee U-'KfU, 

which means that the evolution operator in mass basis can be written as 

U m (L) = e -iL(H m+ A { L)K m ) = fc-iLT^ (15) 

where 4> = e -<i(tr^°m)/3 ) f = J^ m _ (t T Jt m )I/3, and Ji? m = H m + A{L)K m . 

For this case we can thus use the previous expression for T = — (trJ^ n )J/3 by simply replacing A with A(r) 
and then pass to flavor basis by using the U transformation, i.e., T = UTU^ 1 . 

Thus, for any L we can use the spectral decomposition theorem and we find that 

U f {L) = ct>{L)Y, e ~ iLXa{L)p ^ ( 16 ) 
where X a (L) is the ath eigenvalue of T(L) (or T(L)) and 

= ^jttyi m V X ^ L) + Cl(i))/ + K{L)f {L) + f 2(i) l (17) 

3A^(i) + ci(L) L J 

is the projection operator. Everything here is of course L-dependent, since the operator is L-dependent and therefore 
also the eigenvalues. The phase factor is <f>{L) = e -^(tr .^ f (L))/z = e -iL(trje m (L))/3^ 

In the case of a linear matter density of the form A(r) = A + Br, we obtain A(r) = A + Br/2. Similarly, in the 
case of a step function like matter density, relevant to the matter distribution of the Earth, we have 

A 1} < r < Lx, 

A(r) = { A 2 , L x <r<Li + L 2 , , where 2L X + L 2 = L, 

Ax, Lx+L 2 <r< 2Lx + L 2 , 
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which leads to 



A u 0<r<L u 
A(r) ={ A 2 (l-^)+A 1 f 1 L 1 <r<L 1 + L 2 , 

Ax (1 - if) + A 2 if, L 1 +L 2 <r< 2Li + L 2 , 



(18) 



Finally, in the case of an exponentially decreasing matter density A(r) — Ae r / r ° : where A and r$ are parameters 
relevant to the matter distribution of the Sun, we obtain 



A(r) = A— (l - e - r/ro ) 



(19) 



We can see here that in the semi-classical approximation the influence of the density A decays as 1/L with distance 
L from the matter. The evolution should therefore be continued with the vacuum evolution operator as soon as the 
neutrinos leave the matter region. 



IV. DIGRESSION ON THE SEMI-CLASSICAL APPROXIMATION 

Let us consider the semi-classical (s.c.) evolution operator further. It can be written as 

where 

,^\r) =H m r + A^(r)K m . 
Here AW(r) = /J" A(s) ds. Now, the equation of motion for the full evolution operator U m is 



(20) 



i-rU m {r) = Jf m (r)U m (r). 
dr 



This can be integrated to give the equation 



U m (r) =l-i M' m (s)U m {s) ds. 
Jo 

Upon differentiating the semi-classical evolution operator above, we see that, although 

d 



dr 



^\r) = J? m (r) = H m + A(r)K„ 



we can equate i4-U^'{r) with J^C m (r)U^- (r) only when the commutator ,3i?m \r), J^ m (r) 
commutator can be calculated to be 

^«(r),^ m (r)l = / s^(s)ds[H m ,K m }. 

J ,/n ds 



(21) 
(22) 
(23) 

can be neglected. This 
(24) 



Thus, when J Q r s^j(s) ds = ^- (£) for < £ < r is small, the semi-classical approximation to the evolution operator 
is a good approximation to the full evolution operator. For constant matter density this is of course true. For linear 
matter density A(r) = A + Br the coefficient B should be small, i.e., Br <C A or at least Br < A. 
Equation (|22|) can be solved in a systematic way by iteration, leading to 

r-r r-r r-s 

U m {r) = l-i J? m (s) ds + {-if / Jf m (s) / JC(s') ds'ds + ... . (25) 
Jo Jo Jo 

The Hamiltonian Jt? m (r) is the one given in Eq. ( p3| ) . The result to second order is then given by 

2 2 

U m {r) ~ 1 - ir (H m + A{r)K m ) + {-if r - (H m + A(r)K m f + {-if T - (l(r) - A(r)) [H m , K m ] , (26) 
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where 



A(r) = A^(r)/r, A {1 \r) = / A(s)ds, A{r) = 2A (2 \r)/r 2 , A^(r)= A^\s)ds 



By inspection we see that the expression in Eq. (|26|) deviates from an expansion of the semi-classical approximation by 
the terms proportional to the commutator [H m , K m ] and higher order terms in H m and K m . In fact, the commutator 
between the Hamiltonian at different points s and s' is 



[Jf m {s),^ m {s')} = (A(s) - A(s')) [H m ,K m ] = 



dA 
ds 



{s-s')[H m ,K n 



(27) 



which vanishes only for A(s) = A(s'). This is in general true only for constant matter densities, A(s) = A. When 
is large, the contribution of the commutator cannot be neglected. 

We can therefore sum the semi-classical approximation terms and write the solution as 



U m {r) = U s m c -{r) + ai [H m , K m ] + 



(28) 



where 



ai = H) 2 y (A(r)-A(r)) 

and the dots represent higher order terms that vanish when the commutator [H m ,K m ] is neglected. When 
is small, the correction terms are small. 
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V. PROBABILITY AMPLITUDES AND TRANSITION PROBABILITIES 



In the previous sections, we have calculated the evolution operator in the semi-classical approximation. Below we 
will study the corresponding probability amplitudes and transition probabilities. 

The probability amplitude A a p for v a — > vp transition is simply defined as the (/3, a)-matrix element of the evolution 
operator in flavor basis, i.e., 

A afS = (p\U f (L)\a), a,p = e,n,T. (29) 

We now consider transition probabilities for neutrino oscillations in the semi-classical approximation given by Eq. (|l^) . 
Inserting Eq. ( |l6| ) into Eq. (|2|) gives 

3 

A a p = 4>{L) £ e- lL ^P a (L)p a (30) 

o=l 

where 

Pa{L) 0a - (31) 

is the matrix element of the projector P a (L). Here 8 a p is Kronecker's delta. Note that T a p = Tp a and (T 2 ) a p — 
{T 2 )p a - The transition probability P a p for v a — > vp transition is defined as the absolute value squared of the 
probability amplitude A a p. Hence, the transition probabilities in matter are given by the formulas 

3 3 

P a fi = \A a [i\ 2 = S a p - 4 *^2^P a (L) paP^L)^ sin 2 x ab , a, (3 = e, fj,, t, (32) 

a=l b=l 
a<b 

where x ab = (A Q (L) - X b (L))L/2. 
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VI. APPLICATIONS AND DISCUSSION 



The main results of our analysis are given by the evolution operator for the neutrinos when passing through matter 
with varying matter density in the "semi-classical approximation", Eq. (|l6^, and the corresponding expressions for 
the transition amplitudes in Eq. (^0|) and the transition probabilities in Eq. ( 3^) , all expressed as finite sums of simple 
functions in the matrix elements of J^j (or J^ m ) and integrals involving A(r)7the varying matter density. 

As applications, we have calculated the transition probability P^ e for neutrino oscillations for different matter 
density profiles. Our calculations compare two different cases: 

f . An "exact" numerical evolution operator method based on the product of N evolution operators using 

U f (L) = U f (r N - i>r_i)E7>(rtf_i - r w _ 2 ) .. . U f (r 2 - n)t/>(n - r ), (33) 
with the formula (||) used for each step with the appropriate matter density; and 
2. The semi-classical approximation method based on one single evolution operator using 

Uf(L) = 0Ue- iLf U-\ (34) 



In all our examples discussed below, we have used the Earth center crossing neutrino traveling path length, except 
for the last example in which we discuss the Sun. 

Let R® ~ 6371 km be the radius of the Earth and r® ~ 3486 km be the radius of the core with this approximation 
and in a numerical simulation based on Eq. (j|) with a step length of L/N — 2i?®/100 = 127.42 km. The thickness 
of the mantle is then i?® — r® ~ 2885 km with the matter density parameter A\ = A mant i c ~ 1.70 • 10~ 13 eV 
(pi = /Omantic — 4.5 g/cm 3 ), whereas the matter density parameter of the core is A 2 — A corc ~ 4.35 • 10~ 13 eV 

(P2 = Pcore ^ H-5 g/cm 3 ). 

Neutrinos traversing the Earth towards a detector close to the surface of the Earth, pass through the matter of 
varying density densities A(r) where the distances Li, i — 1,2, are functions of the nadir angle h, where h = it — 9 Z ] 6 Z 
being the zenith angle. As h varies from to n/2, the cord L = L(h) of the neutrino passage through Earth becomes 
shorter and shorter. At an angle larger than ho = arcsin(r®/i?®) ~ 33.17°, the distance L 2 = 0, and the neutrinos 
no longer traverse the core. 

The mass squared differences (AM 2 = Am| 2 and Am 2 = Am^) and the vacuum mixing angles (81, 8 2 , 63) used 
here are chosen to correspond to those obtained from analyses of various neutrino oscillation data. We have taken 

AM 2 = 3.2 • 10~ 3 eV 2 , Am 2 = 0, 0i=45°, 6 2 = 5°, 63 = 45°. 

The values of AM 2 and 9\ are governed by atmospheric neutrino data [g7| and the values of Am 2 and 6*3 (LMA) by 
solar neutrino data p8| , where LMA stands for larg e mixing angle (matter) solution. The value of 9 2 is below the 
CHOOZ upper bound, which is sin 2 29 2 — 0.10 |P9| . These choices are the most optimistic ones for obtaining any 
effects in long baseline (LBL) experiments from the sub-leading Am 2 scale [p0| . We should mention though, that 
these data are taken from two neutrino flavor model analyses. 

As a first example, we have investigated the Earth's matter density profile using the published Stacey model for 
the Earth's matter density profile |2(| . The resulting curves of the (exact) numerical evolution operator method for a 
mantle-core-mantle step function approximation of the Earth's matter density profile, the semi-classical approximation 
method, and for reference the (exact) numerical evolution operator method of the Earth's matter density profile are 
shown in Fig. |l|. It has earlier been found by Freund and Ohlsson |Hj that a mantle-core-mantle step function 
approximation of the Earth's matter density profile is a good approximation (even in a three neutrino scenario). 
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FIG. 1. The transition probability P Me as a function of the neutrino energy E v for the mantle-core- mantle step function 
approximation of the Earth's matter density profile. The mean matter density of the mantle and the core were chosen 
to be pmantie = 4.5 g/cm 3 (Anantio ^ 1.7- 10~ 13 eV, Lmantic = 2885 km) and p C orc = 11.5 g/cm 3 (Ax>rc ~ 4.4- 1CT 13 eV, 
icoro = 6972 km), respectively. Parameter values: h = 0, 6 X = 45°, 9 2 = 5°, 9 3 = 45°, Am 2 = 0, and AM 2 = 3.2 • 10" 3 eV 2 . 

The numerical evolution operator method results were carried out using N = 100, i.e., they consist of a product of 
100 evolutions with different constant matter densities for each evolution step, whereas the semi-classical approxima- 
tion method result was obtained with just one single evolution with the average matter density of the Earth's matter 
density profile A§ (p® ~ 7.8 g/cm 3 , which is also the reason why the semi-classical approximation curve only has got 
one resonance peak at E v = E v ,® — 5.4 • 10 9 eV. This peak of course lies inbetween the both resonance peaks (ideally 
at E v = E^ COIC ~ 3.7 ■ 10 9 eV and E v — S v man tie — 9-4 • 10 9 eV) of the exact numerical evolution calculation, since 

^4-mantlc ^4© ^ ^corc- 

In Figs, pi and 0, we have used two different step function matter density profiles. 



0.3 - 



0.2 - 



0.1 - 
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FIG. 2. The transition probability P Me as a function of the neutrino energy E v for a step function matter density profile 
with p x = 4.5 g/cm 3 {A x ~ 1.7 ■ 10~ 13 eV, L x = 2885 km) and p 2 = 5.5 g/cm 3 (A 2 ~ 2.1 • IO -13 eV, L 2 = 6972 km). Parameter 
values: h = 0, 0i = 45°, 9 2 = 5°, 6 3 = 45°, Am 2 = 0, and AM 2 = 3.2 • 1(T 3 eV 2 . 
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Evolution (exact) 

Semi-classit'Lil appmximLiliim 



0.1 - 




io 9 io'° 10" 



FIG. 3. The transition probability P^ e as a function of the neutrino energy E v for a step function matter density profile with 
pi = Ig/cm 3 (Ai ~ 3.8 ■ 10" 14 eV, L x = 2885 km) and p 2 = 2g/cm 3 (A 2 ~ 7.6 ■ 10" 14 eV, L 2 = 6972 km). Parameter values: 
h = 0, 0i = 45°, 02 = 5°, 3 = 45°, Am 2 = 0, and AM 2 = 3.2 ■ I(T 3 eV 2 . 

There appear no double peaks in these figures even though the step function matter density profiles consist of two 
different Af-'s. Furthermore, in Figs. || and]|, the differences between the values of the Ai's and ^'s are the same. 
The step function matter density profile in Fig. |2| could simulate the Earth's matter density profile if the Earth has 
a core, which is much less dense than has been found by geophysics. Note that the absolute error between the two 
curves in Fig. [| is larger than in Fig. |^, whereas the relative error of the curves in Fig. || is larger than in Fig. ||. 

Next, in Figs. H and |l we have studied constant matter density profiles. 



Evolution (exact) 




FIG. 4. The transition probability P^ e as a function of the neutrino energy E v for a constant matter density profile with 
p = 7.8g/cm 3 (A ~ 3.0 • 10~ 13 eV, L = 12742 km). Parameter values: h = 0, 0i = 45°, 6 2 = 5°, 3 = 45°, Am 2 = 0, and 
AM 2 = 3.2 ■ 10" 3 eV 2 . 
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FIG. 5. The transition probability P Me as a function of the neutrino energy E u for a constant matter density profile with 
p = 4.5g/cm 3 (A ~ 1.7 • 10" 13 eV, L = 12742 km). Parameter values: h = 0, 6» x = 45°, 2 = 5°, 6» 3 = 45°, Am 2 = 0, and 
AM 2 = 3.2 ■ lCT 3 eV 2 . 

The constant matter density used in Fig.0 is the average density of the Earth. The curves in this figure could 
be compared with the dotted curve in Fig. [lfln Fig. ||, the constant matter density was chosen to be equal to the 
average density in the mantle of the Earth. In both these figures, the semi-classical approximation method gives an 
excellent agreement with the exact numerical evolution operator method as it should, since in the case of constant 
matter density A = A = const., which means that the two methods arc equivalent. Thus, the very small deviations 
seen in the figures are only due to numerics. 

Then, in Figs. § and ffl, we discuss linear matter density profiles. 




Evolution (exact) 

Semi-chissifLil appmximLiliim 




FIG. 6. The transition probability P Me as a function of the neutrino energy E v for a linear matter density profile with A — 
and BL ~ 3.8 ■ 10" 13 eV (corresponding to p = 5g/cm 3 and L = 12742km if A = -^). Parameter values: h = 0, 9x = 45°, 
82 = 5°, 3 = 45°, Am 2 = 0, and AM 2 = 3.2 ■ 10" 3 eV 2 . 
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FIG. 7. The transition probability P Me as a function of the neutrino energy E v for a linear matter density profile with A — 
and BL ~ 7.6 ■ 10~ 14 eV (corresponding to p = lg/cm 3 and L = 12742km if A = -^). Parameter values: ft = 0, 9x = 45°, 
#2 = 5°, 6> 3 = 45°, Am 2 = 0, and AAf 2 = 3.2 ■ 10" 3 eV 2 . 

The parameter B used in Fig. ^ is larger than that in Fig. ^. The smaller B is the better the agreement between the 
numerical evolution operator method and the semi-classical approximation method becomes. However, linear matter 
density profiles only have theoretical interest, since they are not to be found in Nature at least at large distance 
scales. They could, however, be used on shorter distances scales though, e.g., LBL experiments like K2K, MINOS, 
and CERN-LNGS Jul , where the neutrinos traverse the Earth mantle with cord-like paths. 

Finally, in Fig. ptwe investigated the Sun's matter density profile, which is an exponentially decreasing matter 
density profile. The semi-classical approximation method does not work as well in this case as for step function, 
constant, and linear matter density profiles, since this matter density profile is varying too quickly.P] 




Evolution (exact) 

Semi-dussitLil appmximLilkm 




10 S 



FIG. 8. The transition probability as a function of the neutrino energy E v for the exponentially decreasing matter density 
profile of the Sun with p Q (r) = / 9 (O)e _r/ro , where p (O) = 200g/cm 3 , r = f? Q /10.54 ~ 66000km, and R e ~ 6.96 ■ 10 s m 
IU. Parameter values: 0i = 45°, 8 2 = 5°, 6 3 = 45°, Am 2 = 0, and AAf 2 = 3.2 ■ 10~ 3 eV 2 . 

In conclusion, the semi-classical approximation will be a good approximation for some types of matter density 
profiles. In certain cases, of slowly varying matter density, it is even an excellent approximation. The major advantage 



3 A remark about the neutrino energy interval shown in Fig. |^ is in place. Solar neutrinos have energies of the order of 0.1 
MeV - 10 MeV, i.e., energies much smaller than what is shown in Fig. ^. However, the resonance energy due to the large 
mass squared difference AM 2 = 3.2 ■ 10 _3 eV 2 is given in Fig. ^. Thus, the neutrino energy interval shown in Fig. |§] is only 
of theoretical interest and the figure shows how well the semi-classical approximation works in the region where the transition 
probability changes most. 
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of the semi-classical approximation method as compared to the exact numerical evolution operator method is that we 
only need to calculate one single evolution operator for one single average matter density, the average matter density 
parameter of the considered matter density profile A{L), i.e., we can make the replacement 

71 

U f (L) = l[U f (L i ,A i ) 
v „ ' 

n operators 



iLT 



U f (L)= 

one operator 
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